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Abstract 

Electric vehicles (EVs) offer an attractive long-term solution to reduce the dependence on fossil fuel and 
greenhouse gas emission. However, a fleet of EVs with different EV battery charging rate constraints, that is 
distributed across a smart power grid network requires a coordinated charging schedule to minimize the power 
generation and EV charging costs. In this paper, we study a joint optimal power flow (OPF) and EV charging problem 
that augments the OPF problem with charging EVs over time. While the OPF problem is generally nonconvex and 
nonsmooth, it is shown recently that the OPF problem can be solved optimally for most practical power networks 
using its convex dual problem. Building on this zero duality gap result, we study a nested optimization approach to 
decompose the joint OPF and EV charging problem. We characterize the optimal offline EV charging schedule to 
be a valley-filling profile, which allows us to develop an optimal offline algorithm with computational complexity 
that is significantly lower than centralized interior point solvers. Furthermore, we propose a decentralized online 
algorithm that dynamically tracks the valley-filling profile. Our algorithms are evaluated on the IEEE 14 bus system, 
and the simulations show that the online algorithm performs almost near optimality (< 1% relative difference from 
the offline optimal solution) under different settings. 
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I. Introduction 

Electric vehicles (EVs) are getting more popular as a long-term vehicular technology to reduce the 
dependence on fossil fuel and the emission of greenhouse gases. However, as EVs can consume a relatively 
high amount of power during charging and recharging the batteries, a large fleet of EVs charging at the 
same time may strain the power grid. In fact, with an increase in EV penetration, uncoordinated charging 
can lead to additional power losses and unacceptable voltage variation that overload the power grid [[TJ- 
One way to tackle this problem is to increase the power delivery capacity to deal with new EV peak 
demands, but this will lead to significant infrastructure cost. Another way is to adopt a "smart grid" 
solution, which allows EVs to communicate with the utility that coordinates their charging activities. As 
wireless infrastructure is becoming ubiquitous and more readily available, this solution can be economically 
viable. Furthermore, a coordinated EV charging can improve frequency regulation [2|, smooth out the 
generation intermittency from renewable sources, and increase the efficiency in electricity usage [(3), Q. 

A smart grid is a next-generation power grid that has an overlay of smart meters and communication 
network to enable real-time data exchange between consumers and power providers. As the smart grid 
is purported to achieve better reliability, and greater efficiency in power production and consumption 
[|5J, upgrading conventional grids to smart grids has captured wide attention in the past few years. The 
design of a smart grid thus needs new fundamental theories as well as new technologies in networking 
[|6j. In particular, some residential loads, like washing machines and EVs, might be equipped with smart 
devices that can communicate with the utility company over a network, and dynamically adjust their power 
consumption using feedback control. In this setting, we consider two types of load that are connected to 
the power network: 

• Price-inelastic load: The exact power requested by this type of load must be provided. This corre- 
sponds to standard loads in a conventional grid such as lighting and heating. 

• Price-elastic load: The power delivered to this type of load can vary depending on the current cost 
and a deadline, at which the total energy demand must be satisfied. An example is the charging and 
recharging of EV batteries. 

Considering these two types of loads, the two key problems that we study are: What is the optimal 
charging schedule for EVs to minimize the total power generation cost and EV charging cost? How to 
find a near optimal online algorithm if the future price-inelastic load is uncertain (due to the realistic 
causality constraint)? To formulate these problems, we leverage the well-known optimal power flow (OPF) 
problem and consider its time-dependent extension. 

The solution of the OPF problem optimizes the operation of a power grid, and in general is NP-hard and 
nonconvex [|7J. As the OPF is hard to solve optimally, a number of work are based on local optimization 
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methods, but these often cannot guarantee global optimality (5J. Recently, the authors of [|7), (9j show 
that most practical power configurations surprisingly exhibit a useful property that guarantees zero duality 
gap between the OPF problem and its convex dual relaxation. In particular, [7] suggests solving a convex 
relaxation problem, which is formulated as a semidefinite programming (SDP) problem. When the zero 
duality gap property holds, the relaxation solution coincides with the global optimal solution of the original 
OPF problem. The zero duality gap condition derived in |7J is shown to be satisfied by IEEE systems 
with 14, 30, 57, 118, and 300 buses. This makes efficient polynomial time algorithm to solve the OPF 
problem possible. However, these work look at OPF only as a static optimization problem. 

To incorporate time-varying electricity demand and price-elastic load, we extend the OPF problem to 
a joint OPF-EV charging time-dependent problem that spans over a scheduling period. It consists of a 
finite number of OPF subproblems coupled with one another by the constraints associated with price- 
elastic load, e.g., the EV charging constraints. Similar OPF-related problem formulation has also been 



given in [10]. The joint OPF-EV charging problem is as hard as the OPF problem. In flTOJ, the authors 
prove that the joint OPF-EV charging problem has zero duality gap if the power grid network satisfies 
the zero duality gap condition in [7 ]. The authors propose to solve the convex dual of the joint OPF-EV 
charging optimization, which is a SDP that can be solved optimally in polynomial time. However, this 
approach that uses SDP solvers in [10] is offline, since it assumes that all the future demand for the 
price-inelastic load is known in advance at the start of the scheduling period. Also, for a large number 
of time intervals, this solution is not computationally efficient as the size of a SDP problem grows with 
the number of time intervals. Note that the EV charging problem alone has a scalability barrier, because 
for a long scheduling time, the number of control variables and the possible control sequence can be 
prohibitively large. The dynamic programming approach may not be computationally practical because of 
the huge number of decision states (state explosion). Hence, recent works in EV charging have looked at 



approximation techniques such as approximate dynamic programming in [11 1. 

In this paper, we leverage the zero duality gap result in [|7} to develop both offline and online algorithms 
that solve the joint OPF-EV charging optimization problem. To this end, we propose a nested optimization 
approach that decomposes the joint OPF-EV charging problem into separable subproblems, and then solve 
the decoupled problem using a nonsmooth separable programming technique. We give a novel valley- 
filling profile characterization to the optimal offline solution, and leverage it to design an online algorithm 
for the joint OPF-EV charging problem. The main contribution of the paper are as follows: 

1) For time invariant EV charging cost, we characterize the offline optimal solution to be valley-filling. 
The valley-filling characterization holds true for all network configurations that guarantee the zero 
duality gap condition in the OPF problem. 
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2) With this characterization, we give an optimal EV charging schedule that follows the valley-filling 
characteristic when all the load information are present. This allows us to propose an offline algo- 
rithm that can solve the joint OPF problem with a computational complexity lower than centralized 
interior point solvers. 

3) To account for the causality constraint from the price-inelastic load, we propose an online algo- 
rithm that dynamically tracks this valley-filling characteristic. The online algorithm can be easily 
implemented in a decentralized manner: at each time interval, charging decision is made locally by 
each vehicle after comparing its battery's charging rate limit to an estimated valley level set by the 
utility. We evaluate the performance of our algorithms under different settings of demand and EV 
penetration, and demonstrate that the performance of the online algorithm is near optimal. 

The rest of the paper is organized as follows. The system model is introduced in Section II, and the 
problem formulation is given in Section III. The main analytical results and valley-filling algorithms are 
given in section IV and V. Numerical evaluations are given in Section VI. Finally, the conclusions are 
given in Section VII. We use lower case letters for scalars, blocked lower case letters for vectors, blocked 
upper case letter for matrices, and use (x[t])j to represent the value of the j th component of a vector x[t] 
at iteration t. For clarity of presentation, all the proofs are given in the appendix. 

II. System Model 

We consider a discrete-time model where the time-slot interval matches the timescale at which the 
power grid adjusts its power generation. The goal is to optimize the operation of the power grid over a 
time-interval of interest t £ {0, 1, ... , T}, where T is the scheduling period duration. In practice, T could 
be a day and a slot t could be in the order of minutes. In addition, we assume that price-inelastic loads 
are fixed over each time interval [t,t+ 1]. 

A. Network Structure 

Consider a power network with a set of buses Af :— {1, . . . , n}, a set of generator buses Q C J\f, and a 
set of flow lines C C J\f x J\f. Let z tm and yi„ L be the complex impedance and admittance between bus / 
and m, respectively, and we have yi m = —. Denote Y = (Yi rn , 1 < l,m < n) as the admittance matrix, 

z lm 

where 

yu + Enw yik tfi = m, 

Yim = { -y lm if (I, m) G C, 

otherwise, 
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where yu is the admittance to ground, and n{l) denotes the set of bus k such that (k, I) G C. Now 
let column vectors v = (Vi, V 2 , . . . , Vj^/|) T G C'^' where is the complex voltage at bus k G Af 
and i = (ii, J2, • • ■ , I\Af\) T G C'-^' be the voltage and current vectors, respectively. By Ohm's Law and 
Kirchoff 's Current Law, we have i = Yv. 

From [|7J, [|9j, it is shown that the zero duality gap condition always holds for radial networks, i.e., 
distribution networks. For general networks including transmission networks, the zero duality gap condition 
may hold under mild conditions, e.g., load oversatisfaction [J7J, j9). In addition, this zero duality gap 
condition has been numerically verified for most practical IEEE test cases [J7J, |9). In this paper, we 
assume that the zero duality gap condition is satisfied for a given power grid configuration. 

B. The OPF optimization problem 

The OPF problem finds the optimal operating point that minimizes an appropriate cost function, for 
example generation cost or power loss, subject to certain physical constraints on the power and voltage 
variables p2| . Let us introduce the key notations as follows: 

• Pfc + Qfcj : Apparent power of the load connected to bus k G M (this value is zero whenever bus k is 
not connected to any load). 

• Pk + Qfoj: Output apparent power of generator k G Q. 

• Sim'- Apparent power transferred from bus / G M to the rest of the network through the line (7, m) G C. 

• /fc(Pfc): A convex function modeling the power generation cost for generator k G Q. 

Define (ei, e 2 , . . . , e\M\) as the standard basis vectors in K^L Then the power balance equation is given 

by 

(p* - p{) + (q* - q*)j = v k r k 

= (e^v)(e^i)* = Tr{vv*Y* efce *}. (1) 

By introducing an auxiliary variable W G H^*^ to replace the rank-one matrix vv*, the OPF 
problem can be equivalently formulated as 
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mm V/ fc (p fc ) (2a) 
s.t. Pr n <P fc <^ max , (2b) 

gf n < q fc < Qr x , (2c) 

( ^min )2 < Wfefc < (F max )2; (2d) 

(W a - W^iC < SET, (2e) 

Tr{WY*e fe e*} = p fe - + (q fc - q*)j, (2f) 

W b 0, (2g) 

rank(W) = 1 (2h) 

for every k E Af and (7,m) G £. In the above, the physical limits P fe min , P™ 3 *, Q™™, Qk**> Vj?**, V k max 
and 5™ (where S^™ 3 "* is equal to S™f x ) are given constants. The /*.(•) in the objective function is the 



power generation cost function at bus k. Note that (2b), (2c) are the power generation constraints, (2d) is 



the voltage magnitude constraint, f2e[ ) is the line flow constraint, and pfp is the power balance constraint 
given in ([T]). 

Without loss of generality, we assume that Q = J\f (since when a bus k E J\f is not connected to a 
generator, we let P™ m and p™ ax be both zero). Now, this problem can be relaxed by removing the rank 



constraint ( |2h| ). The relaxed problem without ( |2h| ) is, in fact, an SDP, which is equivalent to the Lagrange 
dual problem of the OPF [|7J. It has been shown in [|7) that, for a large class of practical networks, 
including IEEE 14, 30, 56, 118, 300 bus test systems, the zero duality gap condition is satisfied, and the 
relaxed problem has a rank-one solution. 



C. EV Battery model 

Suppose that each bus k E M can connect to a price-inelastic load and a price-elastic EV battery. 
Furthermore, we assume that each EV battery can absorb or inject only active power at an adjustable 
rate (this can be realized via a smart outlet). For simplicity, we consider that each bus is connected to 
only one EV battery. However, our results can be generalized to the case when multiple EV batteries are 
co-located at the same bus. 

Let (p[£])fc be the control variables corresponding to the active power of the price-elastic load connected 
to bus k E M at time t (these numbers are zero when bus k is not connected to any EV battery). Due to 
the battery specification, an EV at bus k is constrained to charge within a given range [r^r^]. To model 



7 



these charging constraints, we let r k and r k be time-dependent and set (r[t]) k = (r[t]) k = for t that lies 



outside the allowable charging period of an EV at bus k [13|. As such, we have 



(r[t]) k < {p[t]) k < (f[t]) k (3) 

for all k E N and all t E {1, 2, . . . , T}. Let S^, Sfc(O), and rj k denote the battery capacity, initial state of 
charge (SOC), and charging efficiency of an EV connected to bus k, respectively. By the deadline T, the 
EV should be fully charged and this is captured by the total energy stored over the whole time horizon 
Vk (£[*])* At = B k (l- s A (0)). Let c k := B k (l - s k (0))/(r] k At), then the EV charging constraint is 
equivalent to p"3~| : 

T-l 

t=i 

for all k eM. 

Moreover, the battery cannot be discharged when there is no energy available in the battery and cannot 
be charged over its capacity. Hence, for all r E {1, . . . , T — 1}, we have — B k s k (0) < i] k Ylt=i (pM)aA* < 
B k (l — Sfc(O)). By denoting c fc = —B k Sk(0)/(r] k At), the EV charging constraint becomes 



t=i 

for all k E J\f, and r € {1, . . . , T — 1}. Note that if r fc > for all fc, then (|5]) is naturally satisfied. In the 



following, we consider the case when r k > for all k [ 13 1 



III. Problem Formulation 
A. Joint OPF-EV Charging Optimization Problem 

Using the EV battery model and constraints in Sec. |H[ we formulate a joint OPF-EV charging problem 
that optimizes the OPF problem with additional EV charging over discrete time intervals. This is a time- 
dependent OPF charging problem with a coupling of the EV charging constraints over time. Let (p[i])fc 
and (q[t])fc be, respectively, the active and reactive parts of the price-inelastic load connected to bus k E M 
at time t. In the online scenario, we only know the values of p[t] at the beginning of time t. 

Given £ G {1, . . . ,T — 1}, we define p[t] to be the vector of power consumed by the price-inelastic 
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demand, and the joint OPF-EV charging problem can be formulated as 

T-1 T-1 

sjvKn E E /*((pM)*) + E E( a W)*c?w)* ( 6a ) 
{wM,p[ t ]} t=l fc ^ t=l fcgA ^ 

s.t. P fc min < (p[t]) k < Pr x , (6b) 

gr < (q[*])fc < Qr x , (6c) 

(K min )2 < w[f]fefe < (1 / r)2) (6d) 

(W[f] a - WM^JiC < S%*, (6e) 

Tr{W[t]Y*e fc e*} = ((p[f]) fc - (p[t]) fc - (p[t]) fc ) + ((q[t]) fc - (q[t]) fc )j, (6f) 

W[t] h 0, (6g) 

rank(W[t]) = 1, (6h) 

T-1 

EpW = c ' ( 6i ) 

t=l 

r[t) <p[t] <f[t), (6j) 

where («[£])& is a known scalar accounting for the linear charging cost of the EV battery connected to 



bus k at time t, subject to a set of physical limits and power balance constraints. It is shown in [ 10 1 
that, if the duality gap of an OPF problem is zero for a given network, then the duality gap of the joint 
OPF-EV charging problem is also zero in the same network. Hence, solving the relaxed dual problem 
gives the optimal solution to the joint OPF-EV charging problem, provided that p[t] is known for all t. 
However, there will be 0((|jV| + |£|)T) variables in the dual SDP problem, and so the computational 
complexity even when solved by centralized interior point solvers may be prohibitive for large T (the 
computational bottleneck of interior point solvers lies in the inversion of the Hessian matrix when the 
problem size scales up). 

It should be noted that the joint OPF-EV charging problem is coupled through time because of constraint 



(|6Tj). If we can find the optimal p[t], then the problem becomes separable instances of the OPF problem. 
Under the zero duality gap assumption, we will characterize the optimal solution of the joint OPF-EV 
charging problem under the constant price scheme and propose an offline valley-filling algorithm. In 
addition, we propose an online algorithm that solves the joint OPF-EV charging problem by dynamically 
tracking the optimal valley-filling profile. To do this efficiently, we first show how to decouple the OPF 
problem space, i.e., finding W[t], and the EV scheduling problem space, i.e., finding p[t]. 
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B. Decoupling power dispatching from EV scheduling 

While the optimization variables in the Joint OPF-EV charging problem are W[£] and p[t], we see that 
if the charging decision p[t] is known, then all the remaining variables W[i] become separable in time 
t. Hence, solving the Joint OPF-EV charging problem at time t is the same as solving the OPF problem 
with a total demand given by (p[t] + p[t]) as follows: 

F(p[t]+p[t}) :=mm ( £ f k ((p[t]) k ) J (7) 



s.t. (6b),(6c),...,(6h). 



Since the convex dual problem of the OPF can be efficiently solved (7J ||9), we can decouple the power 
dispatching, i.e., finding W[i], from the EV scheduling, i.e., finding p[t], and focus on the following EV 
scheduling problem: 

T-1 

min V F (p[t) + p[t]) + a[t] T p[t] (8a) 

Pit] 

s.t. r[t) < p[t) < f[t) Vt G [l,T- 1], (8b) 

T-1 

J>M=c, (8c) 

t=i 

where F(p[i] +p[i]) in Q returns the optimal value of the OPF problem for a total demand (p[t\ +p[t]). 
Note that a[l], . . . , ot[T — 1] G correspond to the nodal charging cost for each EV at different 
time intervals. In the following, we consider the case when the nodal charging cost is invariant, i.e., 
a[l] = a [2] = . . . = a[T - 1} = a. 

IV. Optimal Offline Algorithm for EV Scheduling Problem 

In this section, we show that F(-) is convex and the optimal solution to EV power scheduling follows 
a valley-filling profile. This characterization allows us to optimally determine the charging rates for EVs 
in 0(1) time. This result therefore reduces the computational complexity of solving ([8]) even when T 
is very large. In addition, we propose an online decentralized algorithm that follows this valley-filling 
characterization and discuss its performance. 

A. Convexity of F(-) 

Recall that the function F(-) given in (|7]) returns the optimal value of the OPF problem for a given 
load demand at each bus. Now, it is not immediately clear whether F(-) is convex or not, because the 
argument to F(-) not only appears in the objective function, but it also appears as part of the constraint 
set. However, by using the zero duality gap condition (when it exists), F(-) can also be expressed as the 
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optimal value of its dual problem. It turns out that for this alternate definition of F(-), the total demand 
(p[t] + p[t\) only appears in the objective function, which leads to the following result. 

Theorem 1. If the zero duality gap condition holds in @, then F : — > R, is a convex function. 

Remark 1. The convexity of F(-) is a direct consequence of the zero duality gap condition of OPF. 

Corollary 1. For the EV scheduling problem in ([8]), since the price-inelastic load (p[l], . . . , p[T — 1]) is 
constant, ((8]) is a convex problem in (p[l], . . . , p[T — 1]), as its objective function is the sum of convex 
functions, and the feasible set is also convex. 

Although F(-) is convex, it can be nonsmooth. In fact, F(-) is a pointwise supremum function in general 
(cf. ( [14] ) in the appendix). In the following, we characterize the optimal solution of ([8]) by making use of 
the convexity of F(-) only, and does not make any assumption on its smoothness. 

B. Characterization of optimal offline solution 

Suppose that the charging rate limits are inactive. Then ([8]) becomes: 

T-l 

mm ^F(p[t] + p[t]) + a T p[t] 
p ^ t=i 

s.t. p[i]+p[2] + ... + p[T-l] =c. (9) 
By convexity of F, we can apply Jensen's inequality and get the following result: 

Lemma 1. If Vi, r[i] = — oo,f[i] = oo, then the EV scheduling problem ((8]) has an optimal solution 
p[l] + p[l] = p[2] + p[2] = . . . p[T - 1] + p[T - 1] = p E + Pl , where p E = (zj=i pW) /(? ~ 1) and 
Pi= (Y,t=xP[t])/(T-l). 

When the charging rate limits are inactive, the optimal solution is a flat profile, i.e., Vi, p[i] + p[t] 
is constant. Next, we consider the case where the charging rate constraints can be active. The optimal 
solution will then no longer be flat, but valley-filling as defined in the following: 

Definition 1. A charging profile is valley-filling, if there exists a unique vector a, such that p[t] = 
[a — p[i]]rw ? Vt, where [x]" = max(Z, min(x,w)). 

In the definition, a can be seen as a valley level that p[t] + p[£] tries to reach unless p[t] is constrained 
by its charging rate limits. A similar definition of valley-filling for EV scheduling can be found in fT4"| . 
However, it should be noted that the valley level in p"4| is the aggregation of all EV demands and is 
a scalar. On the other hand, the valley level in Definition [T] is defined for each individual EV and is a 
vector. Interestingly, the valley-filling characterization is reminiscent of the water-filling notion for power 



allocation to maximize capacity in information theory [15| 
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The next result demonstrates the optimality of a valley-filling profile. It can be proved by using a 
substitution argument, i.e., if there is an optimal charging profile that is not valley-filling, then by convexity 
of F, we can always construct a valley-filling profile that has the same or lower objective value. 

Theorem 2. For a general convex function F(-), a valley-filling profile is optimal to the EV Scheduling 
problem ([8]). 

Corollary 2. A valley-filling profile is a minimizer for any convex function F(-), whether smooth or 
nonsmooth, not just the F(-) defined in (|7]). For example, F(-) can be the I2 norm of the aggregate 
load demand, i.e., F{p\t\ + p[t]) = ($2fceAf((PM)fc + (pM)ai)) 2 - As the total load demand given by 
Ylt=i SfceAf(pM) fc + (pM)fc ?s constant, we can see that the valley -filling profile is also minimizing the 
load demand variance. 

Next, we show the uniqueness of the valley level a. Note that a must satisfy the following for j = 

min{(p[t]),- + (r[t})j} < % < max{(p[t]), + (?[*]) j}, (10a) 

^p[t] = Xi[a-p[t]g]=c. (10b) 

t=i t=i 



If we look at ( |10b| ) component-wise, it is a continuous and strictly increasing function of a^ for slj in 



the box constraint ( |10a[ ). Since ( |10b[ ) is continuous and strictly increasing, we can find a unique a via 
the bisection method for the offline case. This is presented in the following algorithm with e as an error 
tolerance level. 

Algorithm 1 Valley-filling Bisection 

l: Vj,Uj <- max t {(p[t]) j + (f^j}; 

l i ^min i {(p[f]) i + (rM) i }; 
2: for j = 1 -»■ \M\ do 
3: while (||uj — > e) do 

4: m. < _i( u _ ; . + L ; .) ; 

5: if (Eti'K - (p[t])i]g^ > c,) then 

6: Uj ^— rrij; 

7: else 

8: lj rrij; 

9: a «— m. 
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Remark 2. Algorithm [7] determines a in a component-wise manner. Each iteration in the while loop will 
halve the search space for a,, and therefore the computational complexity of the bisection algorithm is 
low. 

Once we have determined the value of a, the scheduling of EV vehicle (solving p[t]) can now be done 
optimally in 0(1) time following the valley-filling characterization in Definition [T] Given the optimal p[t), 
the remaining optimization variables W[£] become separable. The offline algorithm for ([8]) is shown in 
the following. 

Algorithm 2 Offline EV Scheduling 
l: Calculate the valley level a using Algorithm [Tj 

2: for t = 1 ->■ T - 1 do 

3: p[t]Ma-p[t]gj; 

4: Solve the OPF problem with the active load demand set to (p[tj + p[t}); 



Remark 3. Algorithm^decomposes the joint OPF-EV charging problem ^ from a SDP with 0((|jV| + 
|£|)(T — 1)) variables to (T — 1) SDPs each with 0(\Af\ + |£|) variables. Since the complexity of SDP 
interior point algorithms grows superlinearly with respect to the number of variables /[76]/, /[77|/, this 
decomposition leads to a lower computational complexity. 

V. Online Algorithm for EV Scheduling Problem 

Under a causality constraint, we do not assume any knowledge of p[<] until time t. Therefore, we cannot 
use the previous bisection algorithm to find a in an online manner. Instead, we propose an algorithm that 
estimates the valley level, which is denoted by a'[t] and adjusts it dynamically in an online fashion 
(Algorithm [3] below). 

For Algorithm [3j line [T] initializes the first estimation of the (unknown) valley level a. From Lemma 
[T] the ideal valley level is indeed (pe + Pi) if the charging rate constraints are not active. We know 
the value of pe, which is just the charging target c divided by time (T — 1), but the value of pi has 
to be estimated, possibly by learning from historical record of the price-inelastic load. It will be shown 
later that, the more accurate the first estimation is, the better Algorithm [3] performs. Line [3] follows the 
valley-filling characterization outlined in the previous section. However, as the valley level estimation is 
not perfect, we need to take extra steps from line [5] to line [8] to ensure the feasibility of the solution. The 
rationale for line [5] to [6] is to ensure that the charging profile (p[l],...,p[T — 1]) will not overcharge the 
EV batteries at any point in time. Roughly speaking, it means that "if after this instance, even charging at 
the minimal rate will eventually overcharge the EV batteries, then slow down the current charging rate." 
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Algorithm 3 Online EV Scheduling 



1: 


a' [11 +- p E + pi; (pi 


is an estimation of pi) 


2: 


for t = 1 -»■ T - 1 do 




3: 


pW^[a'[t]-p[Cl; 

X |_ J L LJ 1 L JJrl 




4: 


for j = 1 -> |jV| do 




5: 


if (El-i(P^)i > Ci 


" Ef-lli(rM)j) then 


6: 


\V\. L \)3 ~ ^3 Z-ul= 


^~ ,(t\1X\^ — V"^ \ (fi\k]} ,:' (nrevent overcharging 


7: 


if (El-i(pW)i < ^ 


-Ef-m(f [/]),) then 


8: 


(f>[t])j <r- Cj - Y.I 


4+1 ( F W )j _ Efc=\ (pfc] )i ; (prevent undercharging) 


9: 


if (t < T - 1) then 




10: 


(a'[t + l])j <- (a![t 


n T-i-t 


11: 


Solve the OPF problem 


with the active load demand set to (p[i] + p[i]); 



Line [7] to [8] are based on similar rationale to prevent undercharging. Lastly, the estimation of the valley 



level is updated from line M to 10 



Remark 4. Algorithm [5] can be efficiently implemented in a decentralized manner: at the beginning of 
iteration t, the utility broadcasts the estimated valley level a![t], and each EV runs line ^ to line [#] and 



replies with its charging power. The utility collects p[t] and runs line |9j to line 10 to calculate the next 
estimated valley level a![t + 1]. 

The following result demonstrates the feasibility of the output from Algorithm [3} 

Theorem 3. If the charging rate constraints (fx, . . . , f[T — 1]) and (r l5 . . . , r[T — 1]) permit a feasible 
solution, then (p[l], . . . , p[T — 1]) obtained from Algorithm [5] is a feasible solution to ([8]). 

Remark 5. Once the condition in line [5] is triggered, the subsequent charging rate will stay at the 
minimum. Similarly, once the condition in line [7| is triggered, the subsequent charging rate will stay at 
the maximum. 

While line [5] to line [7] ensure feasibility, we may not want to trigger those "if" conditions, because once 
any of them holds true, it means that there is no room left for optimization. Line [9] to line 10 achieve 



this by updating the valley level estimation &'[t] based on information up to iteration t. If the sum of the 
EV load and the price-inelastic load cannot meet the current valley level estimation a'[£], then the next 
estimation &'[t + 1] will be adjusted in the opposite direction to ensure that the "total area" under the 
curve stays the same. This is illustrated in Fig. [2] 

As shown in Fig. [2j the valley level estimation decreases when the sum of the EV load and the price- 
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inelastic load exceeds the current estimation; it increases when the sum is below the estimation, and it 



stays the same when the sum meets the estimated level. This behavior is dictated by line 10, which updates 
the estimation in order to spread out the current error to subsequent estimated valley levels. The result 
below shows that with this dynamic adjustment, the "if" conditions in line [5] and line [7] will be inactive 
in most cases when the first estimation of the valley level is sufficiently good. 

Theorem 4. Assuming that the estimation ofpi is perfect, and p[T — 1] = a'[T — 1] — p[T — 1], then the 
charging profile (p[l], . . . , p[T — 1]) obtained from Algorithm^without Hnel5\tofflis a feasible solution 
to ©. 



Remark 6. Theorem^fjshows that as long as a'[T — 1] — p[T— 1] falls within the range of [r[T— 1], r[T— 1]], 
Algorithm [5] will produce a feasible charging profile without triggering the "if" conditions in line [5] and 
[7| Assuming that the charging rate window [r[T — 1], r[T — 1]] is proportional to the total EV load, this 
suggests that, when there is a high EV penetration level, the dynamic adjustment of the valley level a'[t] 
alone is enough to produce a feasible result. 



Remark 7. As the online adjustment of the valley level (line [9] to line 10) of Algorithm [5] is trying to 



make p[l] + p[l] + . . . + p[T — 1] + p[T — 1] = (T — l)a' 1 , the accuracy in estimating a' x has a direct 
impact on the performance of Algorithm^ 

Fig. [3] compares the different EV charging profiles produced by Algorithm [2] and Algorithm |3j As 
shown, both exhibit the valley-filling characteristic, but the valley level of Algorithm [3] (grey bar) is 
changing dynamically. 

VI. Numerical Results 
Consider the IEEE 14-bus system depicted in Fig. [4j where the circuit specifications and the physical 



limits are given in the library of the toolbox MATPOWER [ 18 1. The system has five generators connected 
to buses 1, 2, 3, 6, and 8. Assume that each of the nongenerator bus 4, 5, 7, 9, 10, 11, 12, 13, and 14 
is connected to an EV load. Enumerate the batteries of these vehicles as 1, 2, . . . , 9. Consider that all the 
batteries are plugged in at time t = 1 and must be fully charged by time t = 25, the charging rate of 
each battery can be controlled only at the discrete time instants 1,2,..., 24. 

Aside from the elastic EV loads, suppose that each bus k £ {1,2, ...,14} is also connected to a 
price-inelastic load given by 

, r l(t) x P k 

(jp[t]) k = W \ t = l,2,...,24, (11) 
l[t) 
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where (Pi, . . . , P u ) is equal to the load profile for a IEEE 14-bus system [19], l(t) follows the average 
residential load in the service area of SCE at different times of the day, (cf. SCE website p0|), and 
l(t) = ^2l =1 l{t). The goal is to optimize the controllable parameters of the power network such as the 
active power supplied by a generator or the charging rate of a battery, which can be modified only at the 
time instants 1, 2, . . . , 24. To this end, we aim to minimize the following cost function 

24 24 

££(pM)* +££«(&[*])*■ < 12 ) 

t=i fceAf t=i keAf 

This cost function has the following features: 

• The generation cost is the total active power generated by all the generators over the time horizon 
from t = 1 to t = 24. 

• The pricing vector of each battery is assumed to be independent of its bus number and invariant over 
time, and we let a = 2 in the following. 

To compare Algorithm [3] with the optimal solution of ( |T2| ) solved by the SDP relaxation approach in 



[10 1, we used CVX pTJ, which is a computational package for solving SDP. Three scenarios will be 



considered in the sequel, and a working example will be shown in Section VI-D 



A. Effect of EV penetration 

In this example, we vary the EV load to be from 10% to 100% of the price-inelastic load, and compute 
the percentage difference. The percentage difference is given by (p* n i ine — Poffline)/Poffline x 100, where 
Poffline an d Poniine are me optimal values obtained by Algorithm [2] and Algorithm [3] respectively. Fig. [5] 
shows the simulation results using three different 24-hour load profiles taken randomly from the SCE 
residential load data [20]. Firstly, we observe that Algorithm [3] solves the joint OPF-EV charging problem 



(|6]) almost optimally. From the three randomly chosen load profile, the worst performance is less than 
0.016% different from p* ffline . Secondly, we can see all three plots go up initially and then eventually 
decrease. This is because at the initial, the EV load is relatively insignificant, and thus Algorithm [2] and 
[3] perform almost the same as there is little to optimize. As the EV load becomes more significant, the 
performance gap grows because Algorithm [3] lacks perfect knowledge. However, a higher EV penetration 
will also lead to a larger room for optimization. Hence, the performance gap decreases and eventually 
approaches zero as the EV penetration increases. 

B. Effect of controllable charging window size 

In this example, we set the charging rate constraint r[t] = 0, V£, hence, the charging window size corre- 
sponds to the upper charging rate limit. We plot the graph for the cases when f[t) equal to 5%, 10%, and 20% 



16 

of the actual EV load respectively. Here, only the residential load profile from 15:00 on Aug. 27th to 
14:00 on Aug. 28th is used. From Fig. [6] we can see that a larger allowable charging rate window leads 
to better results, as the curves corresponding to 10% and 20% charging rate window is below that of 
5% charging rate window. However, once the charging rate window is larger than 10%, any additional 
increase in charging rate window size does not add significantly to performance. 



C. Effect of average load prediction accuracy 

In Fig. |7J we let the initial estimation of the total average load to be a'[l] = (p~E + pi) x (1 + 10z%), i = 
1, . . . , 10, which corresponds to 10%, 20%, . . . , 100% error in estimation of the total average electricity 
load. We set the charging rate window to be 30% of the total EV load. Here, only the residential load 
profile from 15:00 on Aug. 27th to 14:00 on Aug. 28th is used. We can see from Fig. [7] that the more 
inaccurate the initial estimation is, the worse is the performance of Algorithm [3] Also, with a higher EV 
penetration level, Algorithm [3] becomes more sensitive to error in initial estimation. It should also be 
noted that, when the EV penetration level is 80%, an initial estimation error greater than 90% in this case 
causes Algorithm [3] to schedule excessive EV load to be charged at the initial first few time slots that 
leads to the OPF in |7]), i.e., finding W[i], becoming infeasible. 



D. A working example 

In this example, the price-inelastic load variation (see Fig. [8]) is based on the residential load profile 
from 15:00 on Aug. 27th to 14:00 on Aug. 28th (taken from the SCE website (20)). Given the price- 
inelastic load, the EV charging profiles produced by Algorithm [2] and Algorithm [3] are shown in Fig. [9] 
and Fig. [10] respectively. The EV penetration level is set to 50%, and the charging rate window is 20% of 
the EV load, and we assume that there is 10% error in overestimating the initial valley level for Algorithm 
|3j We can make several observations from Fig. [9] and Fig. 10 



1) The charging rate of each EV is high when the inelastic demand (shown in Fig. [8]) is low, which 
corresponds to the valley-filling characterization. The charging profile for EV 3 is flat because there 
is no price-inelastic load present in that bus in the IEEE 14 bus test case specification. 

2) An over-estimation of the initial valley level causes Algorithm [3] to charge the EV batteries faster 
than optimum. As a result, all the EV batteries are fully charged two time slots before the deadline. 

3) The result of Algorithm [3] is still very close to the optimal value, as the percentage difference is 
small: Qo* nline — Poffline)/.Poffline x 100% = 0.0627%. Hence, in terms of power loss minimization, 
Algorithm [3] performs nearly optimally at this setting. 
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E. Runtime comparisons 

In this section, we compare the computational time for the SDP optimization approach that uses SDP 
interior point algorithm to solve the Joint OPF-EV charging problem in [ 10 1 with that of Algorithms [2] 
and [3} The simulation is run on the IEEE 14 bus system, for 6, 12, 24, 48 timeslots respectively^ The 
computational time measured is the average of running the respective algorithm for ten times. 

TABLE I: Performance comparisons 



12 



24 



48 



96 



SDP Optimization 



6.0 s 



13.1 s 



31.5 s 



84.0 s 



262.6 s 



Algorithm 2 



5.8 s 



11.6 s 



22.9 s 



45.8 s 



87.5 s 



Algorithm 3 



5.9 s 



11.6 s 



22.8 s 



45.7 s 



87.4 s 



From Table |I} we can see that the time complexity of Algorithm [2] and Algorithm [3] are comparable. 
Also, both Algorithm [2] and Algorithm [3] have lower time complexity as compared to the SDP optimization 
method in [10|, and the saving in computational time from using Algorithm [2] or Algorithm [3] is more 
significant as T increases. This demonstrates the advantage of the nested optimization approach to solve 
the joint OPF-EV charging problem. 



VII. Conclusion 

We studied a time-dependent OPF charging problem that optimized jointly the operation of the power 
grid and the charging activity of Electric Vehicles. We proved that this problem is convex with respect to 
the total electricity demand, characterized the valley-filling charging profile to be optimal under constant 
electricity price, and proposed a decentralized online algorithm that followed this characterization. At 
each iteration of the online algorithm, each EV calculated its own charging rate according to the valley 
level broadcast by the utility, and the utility guided their charging rate by updating the valley level. The 
online algorithm can be decentralized and thus requires low communication and computation capability. 
Simulation results showed that the online algorithm performed almost optimally in minimizing power loss, 
and the optimal value of the online algorithm approached to that of the offline solution as the penetration 
of EVs increased. However, a higher EV penetration will also lead to a higher sensitivity on the accuracy 
of estimating the average price-inelastic load. 

In this paper, the online algorithm considers a time invariant pricing scheme. That is, the nodal electricity 
price remains constant throughout the scheduling period. However, when there are renewable sources, 

'We used MATLAB version 7.6.0.324 (R2008a) and Windows 7 OS, and ran the programs on an Intel Xeon CPU 2.80GHz machine. 
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electricity prices can vary in real time. In addition, EVs may require charging at different times in a 
more dynamical setting. Incorporating real time pricing, modeling EV arrivals as random events and 
accounting for the additional uncertainties with these extensions are interesting directions for future 
research. Furthermore, our results can be generalized to more complex joint OPF-EV charging problems 
with security constraints, e.g., ramp rate constraints on the power generators formulated along the line in 
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IX. Appendix 
A. Analysis of Valley-filling Characterization 

In this section, we prove Theorem [T| Lemma [T] and Theorem [2] respectively. 

Proof of Theorem 1: Let pE[t] + pi[t] be the argument to F(-) in @. From p\, with the zero duality 
gap condition, -F(p[t] + p[t]) is equal to the optimal value of the following problem: 

maximize h(& , p[t] + p[t]) (13a) 
subject to A(@) h 0, (13b) 

where © is the set of dual variables corresponding to the different physical constraints in ([6]) without 
the rank constraint ( |6h"l ). h(®,p[t) + p[t}) is the Lagrange dual function (in fact, an affine function in 
pM +pW) f° r every feasible ©, and A(&) is a linear matrix inequality constraint in © (cf. the definition 



of Optimization 4 in |7J). Let D denote the feasible set of © that satisfies the constraint ( |13b[ ). Then, 

F(p[t) + p[t)) = su Pg)6D (/i(@, p[t] + p[t])), (14) 



which is a convex function in (p[t] +p[t]) as it is a pointwise supremum of convex functions [23]. □ 
Proof of Lemma^ Ylt^i aT P[A = ot T c, which is a constant. Also, by Jensen's inequality, we have 
F(p[l] + p[l]) . . . + F(p[T - 1] + p[T - 1]) 



>F 



T- 1 

p[l]+p[l] + ... + p[T-l]+p[T-l] 



T-l J \T-1 

where the equality holds when p[l] + p[l] = . . . = p[T - 1] + p[T - 1] = (c + d) /(T - 1). □ 
To prove Theorem [2j we start with the following definition and lemma: 

Definition 2. Denote a > (b,c) > d to mean both a > max(6, c) and d < min(6, c). 

Lemma 2. Let f to be a convex Junction, and a + d = b + c and a > (b, c) > d, then /(a) + f(d) > 
fib) + /(c). 
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Proof: Without loss of generality, assume that a > b > c > d, as a + d = b + c, then there exists a 
k G (0, 1) such that ka + (1 — k)c = b, and kd + (1 — = c. By convexity of /, we have 

kf(a) + (l-k)f(c)>f(b), (15a) 
fc/(d) + (l -*)/(&)> /(c). (15b) 
Summing ( |15a[ ) and ( |15b| ) completes the proof. □ 



Proof of Theorem [2| Consider all the arguments to F(-) to be scalars first, we can get the proof for the 
vector case by repeating the argument for each component of the vector. Suppose (j3'[l], . . . ,p'[T — 1]) is 
an optimal solution to the scalar case of ([8]), and let (p[l], . . . ,p[T — 1]) be the valley-filling load. 

Suppose we have some p'[j] ^ p[j], without loss of generality, assume p'[j] > p[j], as both charging 
profiles are feasible, they need to sum up to the same battery capacity c according to ( [8c] ), hence there 
must exist some k, such that p'[k] < p[k] . From Definition [T] we can observe that the valley-filling 
profile has the minimal deviation from the flat profile, i.e., for ([8]), if (p[l], . . . , p[T — 1]) is a valley- 
filling profile, and (p'[l], . . . , p'[T — 1]) is any feasible charging profile, then a.j — ({p[i))j + (p[i])j) < 
a i - ((p'[«])j + (p'[*])i)> Vi,j. Hence, 

((P\i\)j + (P[«1)i) < (^W)i + (P'WW =► ((P'W), + (P'[«1) 4 -) < a* (16) 

((^w)i + (pM)i) > (m\h + (p'ww =► ((p'W), + (p'M),-) > (i7) 



From ( [16] ) and ( [17] ), we have 

P'\j] + Pbl > a, p'M + < a. (18) 

We consider two cases, p'[j] — p[j] < PE[k] — p'[k] and p'[j) — p[j] > p[k) — p'[k] separately. 
Case 1, suppose that p'[j] — p[j] < p[k) — p'[k], meaning that p'[j] deviates from the valley-filling load 
by a smaller amount. In this case, we swap p'[j] and p'\k\ by p"[j] and p"[k] defined as follows: 

P" [J] = P[j] , f [k] = P' [k] + P' [j] - p[j] ■ 
Now, p"[k] is feasible because p'[k] < p"[k] < p[k] by assumption. Furthermore, p"[j]+p"[k] = p'[j]+p'[k]. 



From ( [16] ), p7] ), and ( [18] ), we have 

p'[k\ +p[k] < (pE[k]+p[k],pE\j\ +p[j\) < p'\j\ +P[j}- 
Hence, by Lemma [2] we have 

F{p[j}"+P[j}) + F(p[k}" + p[k}) < F{p'\j] +p[j]) + F(p'[k]+p[k}), 
which means that applying such swapping procedure will not increase the objective value of ([8]). 



21 



Case 2, suppose that p'[j] — p[j] > p[k] — p'[k]. Then in this case apply the swapping procedure as 
follows: 

P" [j] = P \j) + P i k ] -p[k], P" [k] = Pe [k] ■ 

It can be shown via a symmetric argument to case 1 that, applying such a swapping procedure will not 
decrease the objective value of ([8]) as well. 

Note that each time a swapping is applied, there will be at least one fewer p'[j] that deviates from the 
valley-filling load. Applying a finite number of this swapping procedure completes the proof. 

For the vector case of ([8]), repeat the above argument to each component of the vectors (p[l], . . . , p[T — 
1]). □ 

B. Analysis of Algorithm^ 

In this section, we prove Theorem [3] and Theorem [4] respectively. 

Proof of Theorem l3\ Firstly, the charging profile (p[l], . . . , p[T — 1]) will not overcharge EV batteries, 
i- e -' Sl=i P[k] < c Vt, by line [5] to |6j Secondly, it will not undercharge, i.e., Ylt=i pM — c because of 
line [7] to [8} Hence, Ylt=i pM = c ' me charging sum constraint is satisfied. Thirdly, assuming the charging 
rate constraints permits a feasible solution, then r[t) < p[t) < f[t), Vt. This can be proved by contradiction, 
suppose the charging profile (p[l], . . . , p[T — 1]) obtained from the above algorithm has some such 
that (p[t])j < (r[t])j, note that (p[t])j can only go below (r[t])j if the condition in line [5] is true. Pick 
the smallest t where the charging rate constraint is violated. If t — 1, then that means J2t=i (^]) j > c j 
by line |5} i.e., the charging limit does not permit a feasible solution, which is a contradiction. If t > 1, 
then by line 6 t we have (p[t])j = c,- — J2i=t+i(^\)j ~ Sfe=i(pW)i < (*[t])j , and after rearranging, we 
have X/*=i(p[^])j > c j ~ but that is not possible, because at the end of iteration t — 1, line 

[5] to line [6] will always make sure that ^t=i(p[^])j — c j ~ Sa^CcKDi- The other case where the upper 
charging rate constraint is violated can be proved in a similar way. Therefore, when both the charging 
sum and the charging rate constraints are satisfied, the charging profile produced by Algorithm [3] is a 
feasible solution to ([8]). □ 

Proof of Theorem [?| To prove that (p[l], . . . , p[T — 1]) is a feasible solution to ([8]), we only need to 
show that ESCpW+pM) = (T-l)(p E + Pi) = (T-l)a'[l]. Note that p[T- 1] + p[T- 1] = a'[T- 1]. 
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If T = 2 then we are done, otherwise by line 10 of Algorithm |3j we have 

p[T - 2] + p[T - 2] + p[T - 1] + p[T - 1] = 2a'[T - 2], 
p[T - 3] + p[T - 3] + . . . + p[T - 1] + p[T - 1] = 3a'[T - 3], 



p[l] + p[l] + . . . + p[T - 1] + p[T - 1] = (T - l)a'[l]. 



Hence ^^/(pW + pW) = — l) a '[l]» which is what we set out to prove. 



□ 




8:00 
Line of Day 



12:00 



16:00 



Fig. 1: Base demand curve is the average residential load in the service area of Southern California Edison 
from 20:00 on Feb. 13th to 19:00 on Feb 14th, 2011 p0|. 



23 



Valley Level Estimation 
Ideal Valley Level 
Inelastic Demand 
Total Load 



IT 



IT 



10 15 
Line of Day 



20 



25 



Fig. 2: Illustration of dynamic estimation of valley level. 
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Fig. 3: The base demand curve is the average residential load in the service area of SCE from 20:00 on 
Feb. 13th to 19:00 on Feb 14th, 2011 [20|. The valley-filling curve is obtained using Algorithm [2j The 
valley-filling online curve is obtained using Algorithm [3| 
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EV load as a percentage of the inelastic load(%) 

Fig. 5: Profile 1: residential load from 10:00 on Jul. 6th to 9:00 on Jul. 7th; Profile 2: average residential 
load from 15:00 on Aug. 27th to 14:00 on Aug. 28th; Profile 3: residential load from 1:00 on Mar. 11th 



to 0:00 on Mar. 12th. All load profiles are taken from SCE website [20 1 
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Fig. 7: Effect of percentage error in total average load prediction. 
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Fig. 9: Optimal EV charging rates for the offline case, p* ^ = 79.5348pu. 
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Fig. 10: Optimal EV charging rates for the online case, P* on ^ nQ = 79.5847pu. 



